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ABSTRACT 

We present results of axisymmetric magnetohydrodynamic simulations of the interaction of a rapidly 
rotating, magnetized star with an accretion disk. The disk is considered to have a finite viscosity and 
magnetic diffusivity. The main parameters of the system are the star's angular velocity and magnetic 
moment, and the disk's viscosity and diffusivity. We focus on the "propeller" regime where the inner 
radius of the disk is larger than the corotation radius. Two types of magnetohydrodynamic flows have 
been found as a result of simulations: "weak" and "strong" propellers. The strong propellers are char- 
acterized by a powerful disk wind and a collimated magnetically dominated outflow or jet from the star. 
The weak propellers have only weak outflows. We investigated the time- averaged characteristics of the 
interaction between the main elements of the system, the star, the disk, the wind from the disk, and the 
jet. Rates of exchange of mass and angular momentum between the elements of the system are derived 
as a function of the main parameters. The propeller mechanism may be responsible for the fast spinning- 
down of the classical T Tauri stars in the initial stages of their evolution, and for the spinning-down of 
accreting millisecond pulsars. 



1. INTRODUCTION 

This work studies the interaction of a rapidly rotating 
magnetized star with an accretion disk under conditions 
where the corotation radius of the star, r c = (GM /SI*) 1 / 3 , 
is less than the inner radius of the disk, r^. The radius Td is 
determined by the star's magnetic field and angular veloc- 
ity as well as the mass accretion rate of the disk M and the 
disk viscosity and magnetic diffusivity. Disk accretion is 
disrupted for radii <; Td- At distances ~ Td the disk plasma 
acquires additional angular velocity due to "friction" with 
the magnetosphere of the star which rotates faster than 
the Keplerian motion of the disk. The behavior of the 
matter at Td depends on the ratio between angular veloc- 
ity of the star and Keplerian angular velocity at r^. The 
matter tends to be expelled from the disk if the centrifu- 
gal force is sufficiently larger than the gravitational force. 
This regime of interaction between a magnetized star and 
an accretion disk is called "propeller" regime (Illarionov 
& Sunyaev 1975). This regime was investigated analyti- 
cally (e.g., Davies, Fabian, & Pringle 1979; Lovelace, Ro- 
manova, & Bisnovatyi-Kogan 1999, hereafter - LRBK99; 
Ikhsanov 2002; Rappaport, Fregeau, & Spruit 2004; Eksi, 
Hernquist, & Narayan 2005), and numerically ( Romanova 
et al. 2004b - hereafter RUKL04; Romanova et al. 2005 - 
hereafter RUKL05). 

In the disk-magnetosphcre interaction the magnetic field 
plays the main role in two processes. Firstly, it disrupts 
the inner regions of the disk and thus determines the inner 



radius of the disk. Secondly, the magnetic field determines 
the direction of plasma flow. Matter of the disk at a dis- 
tance r is accelerated in the azimuthal direction if it is 
threaded by field lines connecting it with a star rotating 
more rapidly than the disk; that is r > r c . In the opposite 
case where r < r c , the azimuthal motion of the disk mat- 
ter is slowed down. In the first case the disk plasma gains 
angular momentum from the star and moves outward. In 
second case the disk plasma loses angular momentum and 
moves inward. Thus, in the inner regions of the disk the 
action of magnetic field leads to enhanced accretion rate, 
while in the outer regions to reduced accretion or even to 
outflow of matter. 

The disk-magnctosphere interaction depends essentially 
on the ratio of the inner radius of the disk to the coro- 
tation radius r c . In the model studied here the location 
of the radius rd depends mainly on the magnetic field of 
the star and the disk accretion rate. In the absence of ex- 
ternal plasma the magnetic field of the star is considered 
to be a dipole field with the magnetic moment fj, parallel 
to the rotational axis. Consequently the magnetic stresses 
are proportional to ?'~ 6 . 

For rd <J r c , disruption of the disk is accompanied by 
the formation of "funnel flows" (e.g., Ghosh & Lamb 1978; 
Kbnigl 1991; Shu et al. 1994). Funnel flows were recently 
investigated numerically using 2D and 3D MHD simula- 
tions (Romanova et al. 2002, 2003, 2004a). In this paper 
we consider the opposite limit where rd > r c . In this limit 
one expects that the closed field region of the star will con- 
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tain relatively low-density plasma rotating with the star's 
angular velocity. Outside of this region there is a rela- 
tively dense disk plasma threaded by open magnetic field 
lines. Interaction between the two regions occurs by two 
processes: (1) by matter flow from open field lines to the 
closed field lines as a result of magnetic diffusivity; and (2) 
angular momentum exchange between fast rotating mat- 
ter of magnctosphere and slowly rotating disk. The first 
process can lead to penetration of significant amounts of 
matter with low angular momentum to the closed field 
lines. This can in turn lead to significant deformation of 
the field lines and possibly to their opening. The second 
process may lead to azimuthal acceleration of the disk mat- 
ter sufficient to stop the accretion and to eject matter from 
the disk. An outflow from the system may occur. 

This paper treats disk accretion to a rotating magne- 
tized star using magnetohydrodynamics. For the initial 
set up of the system we divide the space into the "disk" 
and the "corona" with a smooth transition between them. 
The disk is relatively cold and dense while the opposite is 
true for the corona. The disk and corona are separated 
by a specified level of density. In the disk we take into 
account viscosity and magnetic diffusivity of the plasma. 

We consider that both the viscosity and the magnetic 
diffusivity of the disk plasma are due to turbulent fluctu- 
ations of the velocity and magnetic field. We adopt the 
standard hypotheses where the molecular transport coef- 
ficients are substituted with the turbulent coefficients. To 
estimate the value of these coefficients, we use the a-modcl 
of Shakura and Sunyacv (1973) where the coefficient of the 
turbulent kinematic viscosity v t — a v c 2 /£lK, where c s is 
the isothermal sound speed and is the Keplerian angu- 
lar velocity at the given location. Similarly, the coefficient 
of the turbulent magnetic diffusivity rj t — a ( i L c 2 /Q,K- Here, 
a v and ay are dimensionless coefficients which are treated 
as parameters of the model. 

Our numerical simulations have shown that the MHD 
flows which appear as a result of interaction between a 
rapidly rotating magnetized star and accretion disk can 
be divided to two types: "weak" and "strong" propeller. 
In the first case, no outflows are observed (RUKL04), while 
in the second case a large fraction of the mass accretion of 
the disk goes to the outflows (RUKL05). In this paper we 
investigate regime of "strong" propellers in detail. 

The inner disk radius is determined by the disk inter- 
action with the magnetized star. The value r<2 varies with 
time significantly. We can consider however an average 
value of Td- This value depends on the parameters of the 
model which characterize both the star, its mass M, mag- 
netic moment p, and angular velocity f2*; and the disk, 
the coefficients of viscosity and magnetic diffusivity, den- 
sity, temperature and connected with them the accretion 
rate in the disk. The corona may also have a significant 
role. In order to decrease the influence of corona we have 
taken its density to be as low as possible. 

From the parameters of the model we can construct in 
addition to the corotation radius another important length 
termed the nominal Alfven radius, r A = [p 4 / (GMM 2 )] 1 / 7 
(e.g., Davidson & Ostrikcr 1973). Here, M is the disk ac- 
cretion rate at large distances, which can be obtained from 
the simulations. Another characteristic length can be ob- 
tained from the following. Close to the star the magnetic 



field dominates, that is, the magnetic pressure is larger 
than the <fxj> component of the momentum tensor of the 
matter pv 2 +p. Inside the region where the magnetic field 
dominates, the plasma flow is controlled by the magnetic 
field and the flow is along the magnetic field lines. The 
boundary of this region follows from p + pv 2 ^ = B 2 /87r. 
This surface intersects the equatorial region along the cir- 
cle with the radius r m , which we will call the magneto- 
spheric radius. 

The strength of the propeller is determined by the value 
of the dimensionless ratio, rvi/r c or r m /r c , and also by 
other dimensionless parameters of the model such as the a- 
coefficients of the viscosity and magnetic diffusivity. The 
characteristic feature of the disk-magnetosphere interac- 
tion in the strong propeller regime is the formation of an 
outflow of the disk matter along the opened magnetic field 
lines. This outflow is termed the "wind." The matter out- 
flows are strongest in the vicinity of the neutral line of 
the poloidal field, which starts approximately at the inner 
radius of the disk. This is the region of the most intense 
angular momentum transport from the fast rotating mag- 
netosphcre to the disk due to "viscous" stresses. At the 
same time there is a magnetically dominated (Poynting) 
"jet" which flows out along the opened magnetic field lines 
extending outward from the star. The wind has a large in- 
clination angle relative to the rotational axis and for a very 
strong propeller it flows almost along the disk surface. In 
contrast, the jet becomes more collimated as the strength 
of the propeller increases. As the strength increases there 
is a sharp decrease of the accretion rate to the star. 

Section 2 discusses theoretical approach to the problem. 
Section 3 discusses our numerical model. Section 4 dis- 
cusses results and Section 5 summarizes the work. 

2. STRUCTURE AND EVOLUTION OF THE DISK WITH 
MAGNETIC FIELD 

2.1. Non-Magnetic Disk 

For a stationary non-magnetic thin disk, the fluxes of 
matter M and angular momentum L are conserved along 
the disk (see, e.g., Kato, Fukue & Mineshige 1988; Lipunov 
1992) 

M 

rllV r = — — = const , 

L (1) 
r 2 (EV r V<f, + t) = -— = const . 

Here, £ = J dz p is the surface density; t — J dz is 
the "viscous" stress acting in the disk, V r and V<j, are com- 
ponents of the flow velocity. Here, M > is the inward 
mass flux, and L > is the inward angular momentum 
flux. From equations (1) we obtain 

T = ±(Mn-t) = M(n*) . (2) 

2ir y r 2 J 2tt\ r 2 ) W 

Here, SI = V^/r is the angular velocity of the disk, and 
A = L/M > is an angular momentum per unit mass of 
the accretion flow. The viscous stresses are typically very 
small at the inner radius of the disk where the motion de- 
parts significantly from Keplerian. The inner radius of the 
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disk determined by this way is rj s=s X 2 /GM. The "vis- 
cous" stresses are determined by processes of turbulent 
fluctuations which are not fully understood. A simplified 
approach is commonly adopted where the r0-component 
of the viscous stress tensor (integrated over z) is expressed 
as r = — v t TjrdVL/ dr, where v t is coefficient of turbulent 
viscosity. In the region of quasi-Kcplerian disk (O ss Q, K ) 



r = -v t Y>VL K 
2r 



(3) 



E = 



M 



Three-dimensional simulations of the local regions of Kep- 
lerian disks indicate that the viscous stress r constitutes a 
fraction of the gas pressure (integrated over z) in the disk. 
Namely, r = (3a v /2)II = (3a v /2) J dz p, where a v is a 
dimensionlcss coefficient. The value of a v is in the rough 
interval 3 x 1CT 3 - 0.4 (Balbus 2003). Comparing the last 
expression and (3), we obtain the Shakura-Sunyaev formu- 
lae vt = a v Cg/riif, where c 2 . = fl/S is an average isother- 
mal speed in the disk. 

Equations (3) show that the radial accretion speed is 
determined by local properties of the disk. The surface 
density S is set by the local properties as well as the global 
quantity, the accretion rate M. The same is true for the 
nominal inner radius of the disk r j , which is determined by 
the parameter A. To understand how the processes at the 
inner edge of the disk influence the global characteristics 
of accretion process M and L, it is necessary to consider 
non-stationary processes. The continuity equation gives 



as 

~dt 



ld_ 

r dr 



Vsy r ) = . 



Conservation of angular momentum gives 

^ + i|[^ + T)] -0. 
If f2 ~ VIk we obtain 

I d(r 2 T) 



rZVr 



dijdr dr 



2 d(r 2 r) 
dr 



(4) 



(5) 



(6) 



where I = flxr 2 is specific angular momentum. Taking 
the first of equations (3), 



3 v t Q K 

and combining this equation with (4), and (6) we obtain 



8f_ _ 3^t_9 
dt y/r dr 



: df 
dr 



(7) 



where / = r 2 r is the angular momentum flux due to the 
viscous stress. 

The evolution of / = r 2 r is described by a diffusion 
equation with diffusion coefficient ~ v t . Thus, perturba- 
tions generated for example by the star's magnetic field 



at the inner radius of the disk with periodicity At dif- 
fuse out to distances Ar ~ y/ v t At ~ c s \J a v At/flK- If 
Ai<~l/f2/r, then Ar <~ c s ^/o^/Hk < h. Thus, distur- 
bances of the inner part of the disk at say r^ arc not 
expected to influence the more distant parts of the disk 
at r J> 2rd- 

2.2. Analytic Model of Disk Accretion in "Propeller" 

Regime 

Here, we consider the influence on the disk of an aligned 
dipolc magnetic field of a rotating star following the 
treatment of Lovelace, Romanova, and Bisnovatyi-Kogan 
1999(LRBK99). 

The balance of mass, angular momentum and energy 
flows is considered in the region AA'C'C of Figure 1, where 
the surface AA' is along the inner region of the disk, while 
the surface CC is sufficiently far away so that the flow 
can be considered to be undisturbed by the magnetic field 
of the star. The star's mass is M, magnetic moment /x, 
and angular velocity f2» . The accretion disk with accretion 
rate M is located at r > rd and rotates approximately with 
Keplerian velocity. The magnetosphere rotates with angu- 
lar velocity of the star fi» . It is accepted that the angular 
momentum flux through the surface AA' is N — a/j 2 /r 3 , 
where a is a dimcnsionless constant which is an analogue 
of a-coefficient in the model by Shakura-Sunyaev. The 
inner radius of the disk is found to be 



rd J 



7/3 



2a 



7/2 



3 r 2 3/2 



(8) 



Here, r c — (GM/fl 2 ) 1 / 3 is the corotation radius, = 
(li A /GMM 2 ) 1 / 7 is the nominal Alfven radius, k <; I is a 
dimensionlcss constant. (We have simplified the LRBK99 
treatment by taking 6 — in original formulae.) 

Equation (8) was obtained under fairly general assump- 
tions. However, the derivation assumed stationary con- 
ditions and a specific expression for the wind mass, en- 
ergy, and angular momentum outflow was adopted from 
Lovelace, Berk, and Contopoulos (1991). Equation (8) im- 
plies that the inner radius of the disk or effective Alfven 
radius is determined not only by the mass and magnetic 
moment of the star and the accretion rate, but also by the 
star's angular velocity. 

The processes giving the "friction" at the boundary be- 
tween the inner edge of the disk and magnetosphere is 
treated phcnomcnologically by LRBK99 in terms of the 
dimensionlcss coefficient a. The process evidently has a 
turbulent nature and needs a separate detailed investiga- 
tion (see, e.g., Rastatter & Schindler 1999). 

2.3. Disk with Finite Conductivity 

The interaction between the magnetized star and the 
disk depends significantly on the conditions in the corona. 
Magnetic flux of the star extends into the corona and may 
thread the disk. We consider the corona to be a perfect 
conductor. Differential rotation of the foot-points of mag- 
netic field loops threading the star and the disk tend to 
inflate and open unless there are some opposing factors 
which suppress this tendency (e.g., Gold & Hoyle 1960; 
Aly 1980; Lovelace, Romanova & Bisnovatyi-Kogan 1995; 
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Fig. 1. — Geometry of disk accretion to a rapidly rotating star with an aligned dipolc magnetic field. Top panel shows 
radial distribution of angular velocities of the star and of the disk. Bottom panel shows schematic structure of the 
magnetosphere in the propeller regime (from Lovelace et al. 1999). 



Uzdensky 2002). This was also confirmed by a number 
of numerical simulations (e.g., Goodson, Bohm & Winglee 
1999; Fcndt & Elstner 2002). 

One factor is the finite conductivity of the disk plasma 
which leads to reconnection of the inflated field lines (e.g., 
Aly & Kuijpers 1990; Hayashi, Shibata & Matsumoto 
1996; Uzdensky, Konigl & Litwin 2002; Yelenina & Ustyu- 
gova 2005). A second factor is the density of the corona. 
A sufficiently dense corona acts to inhibit opening of the 
magnetic field. Alternatively, a rarefied corona favors the 
opening of field lines with the coronal field passing through 
a sequence of force-free configurations. The finite conduc- 
tivity of the disk leads to the slippage of magnetic field 
lines relative to the disk matter and thus limits the build 
up of the toroidal magnetic field. A dense corona acts to 
oppose twisting of the magnetic field. 

We first consider the influence of the finite conductivity 
of the disk which is assumed to be thin. In a coordinate 
system attached to a point in the disk Ohm's law is satis- 
fied in the form i = AE t , where i is the surface current den- 
sity, E t is the electric field in the plane of the disk, and A is 
the "surface conductivity" (J dza). Thus, in the corona in 
the non- rotating frame, E = — v x B/c, where v is the flow 
velocity in the corona. In the disk, i = A(E + V x B/c), 
where V is the flow velocity in the disk. These relations 
can readily be combined to give 



[(v-V + Cz)xB] t 



(9) 



on the top side of the disk. Here, ( = c 2 /(2nX) is the sur- 
face magnetic diffusivity, z is the unit normal to the disk, 
and index t shows that at the left-hand-side only tangen- 
tial to the disk components of vector are important. If we 
suppose that there is no outflow from the disk (v z = 0), 
then it follows from (9) that a toroidal magnetic field is 



generated at its surface: 



Brt 



Vcf, - V(j> 



(10) 



In order to calculate the toroidal magnetic field we make 
a number of assumptions: First, for the coefficient of mag- 
netic diffusivity we follow the Shakura and Sunyaev a pre- 
scription, rft = (XdC s h, here c s is an isothermal sound speed 
(Bisnovatyi-Kogan & Ruzmaiken 1974). Second, the mag- 
netic field in the disk B z is taken to be the unperturbed 
dipole field of the star, B z = [i/r 3 . Thirdly, the disk is 
assumed to rotate with Keplerian angular velocity while 
the coronal matter above the disk rotates with the angular 
velocity of the star, f2„. Neglecting the variations of the 
tangential electric field E t across the disk, we obtain 



1 _ f h dz _ 
C J-h Vt 



OtdC s 



(11) 



Substituting (11) to (10) we obtain 



Bd 



2r 

CtdCs 



(f2* — Q K )B Z 



2/i 



-(£!* — Hk) 



For ad«l there is a strong twisting of the magnetic field 
linking the star and the disk, 



B,p 

~bZ 



grja. - n K \ 

adC s 



> 1 



if the considered distance r is far from the corotation ra- 
dius where fi* = Qk ■ This equation was derived earlier 
by Lovelace et al. (1991). 

The discontinuity of the azimuthal component of the 
magnetic field at the disk surface z — means that there 
is a radial surface current flowing along the disk. The in- 
teraction of this current with the vertical magnetic field, 
B z , leads to the loss of angular momentum from both sides 



5 



of the disk at the rate —2rB t j ) B z /4n per unit surface area. 
For stationary conditions, the conservation of the angular 
momentum balance gives 



TTa d c s r' 1 



-(ft* — Qk) ■ 



Assuming ad, c s = const, and V$ = Vk{t), this equation 
can be integrated to give 

M( \\ m 2 (2VGM ft*\ 



where A = L/M is a constant of integration. 

We identify the point r = r d , where r = 0, with the 
inner radius of the disk. At this distance the viscous and 
magnetic stresses lead to significant deviation of rotational 
velocity from Keplerian. For a sufficiently large magnetic 
moment and a rapidly rotating star, the deviation from 
Keplerian motion is determined by the magnetic field. We 
can then estimate r d by comparing the first and the last 
terms on the right-hand side (12). This gives 



5/2 



7/2 



MVGMa d c s adh* 



where ta = ((J, 4 /GMM 2 ) 1 / 7 is the nominal Alfven radius 
and /i* = c s /ft* is the hydrostatic thickness of the disk at 
the corotation radius. Thus, r d = r^Tvi/ay/i*) 2 / 5 S> va- 
This derivation suggests that r d ^> r c and r d 3> rt = 
\ 2 /GM. 

In order to estimate the influence of the magnetic field 
on the disk, we consider a tube of field lines with one foot 
point frozen to the star and rotating with the angular ve- 
locity ft*, and the other foot-point in the disk rotating 
with angular velocity ft ^ ft*. The field lines will be de- 
formed due to the differential rotation and this gives rise 
to an azimuthal magnetic field component. Furthermore, 
the field lines tend to open as a result of the differential 
rotation. At the same time the azimuthal component of 
the field acts in combination with the poloidal field applies 
a torque to the disk matter. To estimate the net effect on 
the disk we assume that the opening of the field lines oc- 
curs on a time scale At ~ 1/Aft, where Aft = |ft* — ft|. 
During this time the azimuthal magnetic field builds up 
to a value ~ B z at the surface of the disk. The mag- 
netic stress on the disk during this process exerts a torque 
on the disk (from both sides). The corresponding loss of 
angular momentum (per unit area) from both sides of the 
disk is At rB c f>B z /2ir ~ rB 2 / '(27rAft). Of course angular 
momentum per unit area of the disk is simply r 2 £ft. The 
magnetic field will have a significant influence on the disk 
when 

rB 2 

27rAft ^ 

Assuming the disk to be roughly Keplerian, ft w GM/r 3 

with £ = M/3nv t7 and the magnetic field B z = fi/r 3 , we 
obtain the condition 



ft|ft* - n\ <, 



M 2 /r 



= 3v t H 
2wM/{3TTU t ) ~ 2 Mr 7 



(13) 



for the magnetic field to influence the disk motion. 

In the limits of slowly and rapidly rotating stars, in- 
equality (13) provides an estimate of the inner radius of 
the disk: For a slowly rotating star, (ft* <C ft), 
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Y' 7 = r J£?L\ 



\G 3 M 3 M 2 



1/7 



T A 



\ G 2 M 2 



For a rapidly rotating star, (ft* » ft), 



rd 



/ 2 2 

/ a v c /J, 



\GMAin* 
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Here, we have assumed v t = a v c 2 /ftif, a v = const and 
dropped insignificant numerical coefficients. The first case 
corresponds to the situation where the corotation radius 
r c is much larger than the inner disk radius. The sec- 
ond case corresponds to the opposite inequality. In the 
first case the disruption of the disk occurs in the Funnel 
Flow regime. Matter loses its angular momentum because 
of "friction" with the slowly rotating magnetosphere and 
accretes to the star along the dipole field lines, forming 
the funnel flows. In the second case the disruption of the 
disk occurs in the Propeller Regime where matter acquires 
angular momentum "friction" with the rotating magne- 
tosphere and is driven away from the star, forming the 
outflows (jets and winds). Sketches of the flows are shown 
in the figure 2. Sketches of the radial variation of the disk 
angular velocity in the equatorial plane in the two cases 
are shown in Figure 3. 

The transition at ~ r d between the magnetosphere re- 
gion which rotates with the angular velocity of the star 
ft* and the disk which rotates with Keplerian angular ve- 
locity VIk (more exactly a part of the disk which rotates 
with quasi-Keplerian velocity) is not sharp, but instead 
occurs over some distance where instabilities may have a 
significant role. Instabilities can act generate turbulence 
at the inner edge of the disk. We expect turbulence will 
be strongly developed in the "propeller" regime. Figure 4 
shows the dependence of the specific angular momentum 
I = fir 2 on radius in the two cases. 

One can see that in the "propeller" regime there is a 
region where dljdr < 0, which is the Raylcigh criterion 
for axisymmetric instability (Chandrasckhar 1961). The 
decrease of the angular momentum with the distance is 
only one of the reasons which may give rise to turbulence 
in the disk. Nevertheless, it is reasonable to expect that 
the turbulent viscosity and magnetic diffusivity are signif- 
icantly enhanced at the inner edge of the disk compared 
with the regions of the Keplerian disk at larger distances. 

3. FORMULATION OF THE PROBLEM 

3.1. Physics Problem 

As discussed the nature of the interaction of the magne- 
tized star with the accretion disk depends on many factors. 
Nevertheless, we can try to classify the resulting configu- 
rations basing on the relations between three lengths, the 
corotation radius r c , the magnetospheric radius r m , and 
the inner radius of the disk r d . We recall that for r m we 
have in mind the location in the equatorial plane where 



B 2 

8?r 



(14) 
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Fig. 3. — Sketches of the radial distribution of the angular velocity in the non-propeller (left panel) and the propeller 
(right panel) cases. 




Fig. 4. — Sketches of the radial distribution of the specific angular momentum in the non-propeller (left panel) and 
propeller (right panel) cases. 



It is convenient to use the radius r m for comparison with 
the other two radii even though it depends on the re- 
sults of simulations. We also note that in the subsequent 
parts of the paper we will consider only the "propeller" 
regime. That is, the star is rapidly rotating in the sense 
that r c < Td- Classification based on the ratio between 
r<j and r m is not similarly clear, because all variables in 
the (14) are rapidly varying functions of the position of 
the point. This is why significant variation of the param- 
eters of the system, for example, the magnetic moment of 
the star, may not lead to significant variation of r m . An 
important quantity is sound speed in the corona, which de- 
termines variations of the density and pressure (see section 
on Initial conditions). Results of our numerical simulations 
show that it is always a true that w r m , but in spite of 
that the matter flow may have different morphology. In a 
number of cases the magnetic field forms configuration in 
which some of the poloidal field lines are well-collimated 
along the symmetry axis, while at larger distances the field 



part is pushed toward the equatorial plane (RUKL04) . In 
other cases the collimation along the open field lines is 
less prominent, but matter flows from the disk along the 
neutral layer where B p w (see also RUKL05). We trace 
such qualitative differences in the final picture of the disk- 
star interaction to the initial temperature or sound speed 
in the corona which determines the initial differences of 
density and pressure. 

In the considered model we assume a cold quasi- 
Keplerian disk surrounded by a relatively hot corona. The 
disk and corona are in the gravitational field of the central 
object, a star of mass M which rotates with angular ve- 
locity SI* around the same axis as the disk. The star has 
a magnetic moment fi which is aligned with the rotational 
axis of the star. 

Having in mind application of the results of the numeri- 
cal modeling to T Tauri type stars, we mention typical pa- 
rameters, the mass of the star, M = 0.8M Q = 1.6 x 10 33 g, 
its radius r* = 2r Q = 1.4 x 10 11 cm, the inner radius of 
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the disk, rj = 5r* = 7 x 10 11 cm, the magnetic field at the 
surface of the star (in the equatorial plane), £?* = 10 3 G, 
the temperature of the corona T cor — 10 6 K, and the tem- 
perature of the disk T disk = 5x 10 3 -ftT . 

The Kcplerian velocity at the inner radius of the disk is 
Vk = y/GM/rd w 1.2 x 10 7 cm/s. The isothermal sound 
speed in the corona is c cor — yjRT cor /m w 1.3 x 10 7 cm/s 
(for fully ionized hydrogen m = 1/2); in the disk Cdisk = 
y/RTdisk = 6.5 x 10 5 cm/s, where R = 8.31 x 10 7 ergs/(mol 
K) is the gaseous constant. The magnetic moment of the 
star is fi = B»rJ = 3.43 x 10 38 G cm 3 . We treat the den- 
sity in both the disk and the corona as free parameters. 
For a fixed sound speed in the corona (or the disk) and a 
given magnetic field of the star, the choice of say the disk 
density determines the ratio of the Alfven speed and the 
sound speed, or equivalently the ratio of the magnetic and 
gas pressures. We chose the density such that in at least 
part of corona near the star, the magnetic energy density 
dominates, B 2 /87r > pv 2 ^ + p. At the inner edge of the 
disk the Keplerian velocity is related to the sound speed 
as vk ■ c cor ■ Cdisk "1:1: 0.05. 



ordinates equations (15) have the form: 

dp 1 , 2 x 1 8 , . n s 
ot r z or ' r sin 9 d9 ' 



d(pv r ) | i d , 2t . 

dt r 2 dr " 



dt 



1 o> 
r sin 89 



sm 



Tee + Xt>, 



1 d , 3 _ x . 1 d 



6T r6 ) 



P9r 



djpve) 
dt 



1 d 



1 d 



= cotan6* 



rd6' 



T~d><t> — %e 



= 0, 



- — (r 3 T r4> ) + —^- — ( sin 2 9T e<p ) = 0, 
r i ar \ f/ rsm 9 o0 ' 



3.2. Basic Equations 

We assume that the plasma flows are described by the 
equations of magnctohydrodynamics. Further, we consider 
that any shock waves in the flow are outside of the con- 
sidered region, for example, inside the considered region, 
on the surface of the star. This allows us to use an energy 
conservation equation in the form of an entropy continu- 
ity equation. The system of equations in a non-rotating 
reference frame is: 



dpv 
~~d7 

dB 

!h + 

d(pS) 



dt 



+ V T 
cV x E 
+ V • (pSw) 



, 
P g 
o , 
Q ■ 



(15) 



Here, p is the density and S is the specific entropy; v is 
the flow velocity; B the magnetic field; T the momentum 
flux-density tensor; E the electric field; Q is the rate of 
production or loss of entropy per unit volume; c is speed 
of light; and g = — (GM/r 2 )r is the gravitational acceler- 
ation due to the star of mass M. The total mass of the 
disk is negligible compared to M. 

We consider axisymmetric MHD-flows and use spher- 
ical coordinates (r,0,(j>). The polar angle counted from 
symmetry axis, all derivatives d/dcf) = 0. In spherical co- 
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r dr 



dt 



(rcE e ) 



sm 
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dB 1+ ld_, 
dt r dr 



cE r ) = , 



d(pS) + 1 d 



dt ' r 2 dr 
The stress tensor T is 



(pSv r ) + 



1 d 
r sin 9 d9 



I B 2 B,lh. 
Mk = PViVk + POik + I ~^ d ik j^— 



(sin9pSve) = Q 



Tik — T ik + Tjfe 



Here, p is the pressure; is the "viscous" stress caused 
by the turbulent fluctuations of the velocity and magnetic 
field; g r = —GM/r 2 is the acceleration due to gravity. 
The plasma is considered to be an ideal gas with adiabatic 
index 7 = 5/3, so that S = \n(p/ p 1 ). 

We assume that the stress due to fluctuations of the 
velocity and the magnetic field can be represented in the 
same way as collisional viscosity by substitution of the 
turbulent viscosity coefficient. Moreover, we consider that 
the viscous stress is determined mainly by the gradient of 
the angular velocity because the azimuthal velocity is the 
dominant component in the disk. The dominant compo- 
nents of the tensor T ik in spherical coordinates are: 

^4, = -v t prsm.9— , 



T64> 



^ipsmt 



jdu> 
~d9 ' 



Here u — v^/r sin 9 is angular velocity of the plasma and 
v t is the coefficient of the turbulent viscosity. 
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Separating out the viscous stress gives 

1 9(sin 2 6T e0 ) 



dt 



dr 



13/ 4 • 

r J or \ or J rsm 



r sin U 

1 a 



2 961 



96» 

• 3 



(16) 

where T rr j, and Xg^ are components of the inviscid part of 
the stress tensor. 

The viscosity leads of course to dissipation of the ki- 
netic energy and its conversion to thermal energy and to a 
corresponding increase of the entropy. The rate of energy 
dissipation per unit volume is 



'ik 

2pv t 



— — sin 



dui 
dr 



du> 
~d6 



We assume that the finite conductivity of plasma as well 
as "viscous" stresses are connected mainly with the turbu- 
lent fluctuations of the velocity and magnetic field. The 
induction equation averaged over the small-scale fluctua- 
tions has the form 



9B 

— V x (v x B) + cV x E f 

dt 



. 



(17) 



Here, v and B are the averaged velocity and magnetic 
fields, and = — ((v' x B')) /c is electromotive force 
connected with the fluctuating fields. Because the tur- 
bulent electromotive force E^ is connected with the small- 
scale fluctuations, it is reasonable to suppose that it has 
a simple relation to the ordered magnetic field B. If we 
neglect the magnetic dynamo a-effect (Moffat 1978), then 
((V x B')) = — rftV x B, where r/ t is the coefficient of tur- 
bulent magnetic diffusivity. Equation (17) now takes the 
form 



<9B 
~dt 



V x (v x B) + V x (f/tV x B) = 



(18) 



We should note that the term for E^ formally coincides 
with the Ohm's law 



J = —V x B = w.t 

47T 4:TTr]t 



E T . 



The coefficient of turbulent electric conductivity a = 
c 2 /47r?7 t . The rate of dissipation of magnetic energy per 
unit volume is 

^ = f(VxB) 2 . 
a 47r 

To calculate the evolution of the poloidal magnetic field 
it is useful to calculate the ^-component of the vector- 
potential A. Owing to the assumed axisymmetry, 



B r — 



1 9(sin6»A ) 



1 d{rA^) 



r sin 6 d6 r dr 

Substituting B = V x A to the induction equation gives 
<9A 



dt 



+ cE = V X , 



(19) 



where \ is an arbitrary function. Thus, 

A -(vxB) 4> + Vt (VxB)^ = 0, 



dt 



dA d 



dt * \ r dr 2 



ld 2 (rA <j) ) 1 d 1 9(sin6»A ) 



(20) 
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r* OU sm ( 

The azimuthal component of the induction equation gives 



dt r dr \ dr J r 2 d9 \ sin 9 d6 
_ 1 fd[r(y x B) e ] 9(v x B)/ 



r \ dr d9 
The Joule heating rate per unit volume is 



(21) 
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f d(rBg) dB r \ 2 | f d{rB^) \ 2 | 1 f d(sm6B^ ~ 



V dr 



de 



V dr J sin 2 6< 



3.3. Initial Data 



do 



A variety of simulation runs were calculated for a range 
of values of the star's magnetic moment and its angular 
velocity and for different values of the viscosity and mag- 
netic diffusivity. Changes in the magnetic moment of the 
star (at t = 0) of course changes the location of the surface 
where p + pv 2 = B 2 /87r. 

The initial conditions have a significant role in giving a 
"smooth" transition of the system to possible stationary, 
quasi-stationary or a quasi-periodic regime. This role of 
the initial conditions provides a measure of their validity: 
the time-averaged location of the Alfven surface (in the 
equatorial plane) should not be appreciably different from 
the its initial value. In the opposite significant part 

of the computational time can be spent for establishment 
of more or less equilibrium distribution of the matter and 
the field. 

The initial magnetic field is determined by the magnetic 
moment of the star. The field is a dipole with field lines 
connecting the star and the disk. It is clear that the dif- 
ferential rotation of matter along the field lines may lead 
to the generation of an azimuthal magnetic field and to 
the opening of some of the field lines (e.g., Lovelace et al. 
1995). The qualitative picture of the poloidal field lines is 
shown in Figure 5. The dashed line shows the line along 
which the r-component of magnetic field is zero. A final 
field configuration could be determined and used as an ini- 
tial condition. Nevertheless, we have found it efficacious 
to start with the poloidal dipole field, and let the system 
evolution open the field lines (see also Romanova et al. 
2002). 

Initially the matter of the disk and corona are assumed 
to be in mechanical equilibrium, and the magnetic field is 
the current-free dipole field of the star B r = 2/i cos 6>/r 3 
and Bg = psin9/r 3 . In addition, the initial density dis- 
tribution is taken to be barotropic with 



Pip) 



P IT disk , P> Pb and rsin6» > r b 
p/T cor , p < Pb or rsin6» < r b . 



9 




Fig. 5.— Sketches of the initial (left) and final (right) configurations of the poloidal magnetic field which results from 
the difference in angular velocities of the star and the disk. 



The level surface of pressure pb separates the cold matter 
of the disk from the hot matter of corona. At this sur- 
face the density has a discontinuity from value pb /T disk to 
value Pb/T cor . 

Because the density distribution is barotropic, the an- 
gular velocity is a constant on coaxial cylinderical surfaces 
about the z— axis. Consequently, the pressure may be de- 
termined from 

F(p) + $ + $ c = E = const . 

Here, $ = — GM/\r\ is gravitational potential, <I> C = 
Irsine^ 2 ^^^ 1S centrifugal potential, which depends 
only on cylindrical radius rsin#, and 

{Tdisk ln(p/p b ) , p > p b and r sin6> > r b , 
T cor \n(p/p b ) , p < pb or rsinfl < r b ■ 

Initially, the inner edge of the disk is located at r b = 5 
in the equatorial plane. The angular velocity of the disk 
is slightly sub-Keplerian 0,(6 = tt/2) = (1. - 0.003)fi K , so 
that the density and pressure decrease towards periphery. 
Inside the cylinder r < rb the matter rotates rigidly with 
the angular velocity 0(r 6 ) = (1. - 0.003) y/GM/rf. To in- 
sure the smooth start-up, the angular velocity of the star 
is changed gradually from its initial value Vl cor = 5~ 3//2 ~ 
0.09 to a final value f2» during three Keplerian rotation 
periods at r = 1. 

Below we show results of our simulations for the fol- 
lowing set of dimensionless parameters: the gravitational 
potential is such that GM — 1; the dipole moment of 
the star is p = 10; the angular velocity of the star is 
O* = 1; the initial "temperature" of the disk is Tdisk = 
(p/p)disk = 0.0005; the initial temperature of the corona 
is T cor = (p/p) cor = 0.5. We note that the dimensional 
temperature is T(K) = mv^T/R, where vo is the velocity 
scale, m is the number of barions per particle (including 
free electrons), which by itself depends on the tempera- 
ture, chemical composition and other factors. 

3.4. Transition to the Dimensional Variables 

Results obtained in dimensionless form may be applied 
to objects with widely different scales. Below we give ex- 
amples for classical T Tauri stars (CTTSs) and for neutron 
stars with a weak magnetic field (millisecond pulsars). 



For CTTSs we adopt the following scales. The mass of 
the star is M = M* = O.8M = 1.6 x 10 33 g. The length 
scale is r = 2r st = 4r Q = 2.8 x 10 11 cm. The magnetic 
field at the surface of the star is -B* = 10 3 G. From these 
numbers we obtain the other scaling variables: The time- 
scale to = 0.166 days; velocity vo — 1.95 x 10 7 cm/s; den- 
sity po = 4.1 x 10~ 13 g/cm 3 ; pressure po = 156 erg/cm 3 ; 
accretion rate M = 6.27 x 10 17 g/s = 9.84 x 10~ 9 M Q /yr; 
angular momentum flux Lo = 3.43 x 10 36 gcm 2 /s 2 . For 
these parameters, the magnetic moment of the star is 
p* = 2.74 x 10 36 Gem 3 ; the star's rotation period is 
-P* = 1.04 days, the initial temperature in the disk is Td = 
2,300 K, and coronal temperature is T c = 1.15 x 10 6 K. 
Subsequently we give time in periods of rotation of the 
star. 

For a neutron star, we give an illustration of the scales. 
The star's mass is M = 1.5M = 3.0 x 10 33 g. The 
length scale is r = 2r» = 2 x 10 6 cm. The magnetic field 
at the surface of the star is = 10 9 G. Other scaling 
variables include the time-scale to = 2 x 10~ 4 s, the veloc- 
ity v — 10 10 cm/s, the density p = 1.56 x 10~ 6 g/cm 3 , 
the pressure po — 1.56 x 10 14 erg/cm 3 , the accretion rate 
Mo = 6.25 x 10 16 g/s = 10 -9 M©/yr, and the angular mo- 
mentum flux Lo = 1-25 x 10 33 gcm 2 /s 2 . For these parame- 
ters, the magnetic moment of the star is /i* = 10 27 Gem 3 , 
the rotation period of the star is P* = 1.26 x 10~ 3 s, the 
initial temperature in the disk is Td = 3 x 10 8 K, and the 
coronal temperature is T c = 3 x 10 11 K. 

3.5. Numerical Method 

For the numerical integration of the MHD equations 
including magnetic diffusivity and viscosity, we used a 
method of splitting of the different physical processes. Our 
simulation algorithm has a number of blocks: (1) a hy- 
drodynamic block in which we calculated the dynamics of 
the plasma and magnetic field with dissipative processes 
switched off; (2) a blocks for the diffusion of the poloidal 
and azimuthal components of the magnetic field calculated 
for frozen values of the plasma velocity and thermody- 
namic parameters (density and pressure); and a block for 
the calculation of viscous dissipation in which we took into 
account only r</> and 9(j) components of the viscous stress 
tensor. 

The system of MHD equations (15) was integrated nu- 
merically in the region r int < r < r ext , < 9 < n/2. We 
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used grids which were equally spaced in 0. The steps in 
the radial direction were chosen so that the poloidal plane 
cells were curvilinear rectangles with approximately equal 
sides. A typical region was 0.5 < r < 50, < 9 < n/2 
with the grid 85 x 31 cells. Figure 6 shows the grid in the 
entire simulation region (right panel) and in the inner part 
of the region (left panel). 

In the hydrodynamic block, the ideal MHD equations 
are integrated using an explicit Godunov-type numerical 
scheme. For calculation of fluxes between the cells we 
used approximate solution of the Riemann problem anal- 
ogous to one described by Brio & Wu (1988). To guar- 
antee the absence of the magnetic charge, we calculated 
at each time-step the ^-component of the vector-potential 
Afj,. This was then used to obtain the poloidal compo- 
nents of the magnetic field (B r , Bg) in divergence- free form 
(Toth 2000). 

In the block where the diffusion of the poloidal magnetic 
field is calculated, we numerically integrated the equation 
(20) for the ^-component of vector-potential. During this 
calculation we frozen the values of on the inner and 
outer boundaries of the simulation region. In the equato- 
rial plane we have conditions of symmetry, dA^/dO = 0. 
On the symmetry axis we have A^ = 0. 

Equation (20) was approximated with an implicit dif- 
ference scheme. The approximation was chosen so that 
the operator on the implicit time-layer was symmetric and 
positive. 

For solving the system of equations on the implicit time- 
layer, we used ICCG method (Incomplete Cholesky Con- 
jugate Gradient method). Because the size of the grid 
cells and values of the coefficient of the magnetic diffusiv- 
ity vary strongly in space, the elements of the matrix of 
the system also vary strongly. To remove this undesirable 
property, we changed the matrix so that it has diagonal 
elements equal to unity. 

In the block where the diffusion of the azimuthal compo- 
nent of the magnetic field was calculated, we numerically 
integrated equation (21). At the inner and outer bound- 
aries the B<p was frozen at this computational block. Along 
the rotation axis and on the equatorial plane B^ = 0. 

Equation (21) was approximated by a numerical scheme 
with symmetric and positive operator on the implicit time 
layer. The corresponding system of linear equation was 
solved by the ICCG-method with pre-conditioning. 

In the block of our code where the viscous stress is cal- 
culated, we numerically integrated the equation (16) for 
the angular velocity of matter u> = v^/rshiO. At the in- 
ner boundary of the simulation region we take u> = fi*, the 
angular rotation velocity of the star. At the outer bound- 
ary, lu is taken to be fixed and equal to the corresponding 
Keplerian value. On the axis and in the equatorial plane 
we have the condition of zero stress for the 9(f> -component 
of the viscous stress tensor. 

Equation (16) was approximated by a numerical scheme 
with the symmetric and positive operator on the implicit 
time-layer. The corresponding system of linear equations 
was solved by the ICCG-method with pre-conditioning. 



a v and ay- In the following we discuss results for the 
reference case p, = 10, f2* = 1, a v = 0.2, and ay = 0.2. 

The reference case was calculated out to a time corre- 
sponding to 2500 rotation periods of the star. We observed 
that the process of interaction between the magnetosphcrc 
and the disk is accompanied with (1) oscillations of the 
disk, (2) quasi-periodic opening and closing of the mag- 
netic field, (3) quasi-periodic outbursts to winds and jets. 
This evolution arises from the diffusive mixing of the disk 
matter with the magnetic field of the star, the opening of 
the magnetic field lines linking the star and disk, unload- 
ing of the disk matter to the star and to the outflows, and 
the outflow of angular momentum from the star and from 
the disk. 

The left-hand panel of Figure 7 shows the variation of 
the inner radius of the disk which was calculated at the 
density level, po = 0.3. Note that near the disk-corona 
boundary, p(r, 9, t) = p , the density gradient is very large 
so that the choice of po is not important. During an ini- 
tial time t ~ 100, part of the magnetic flux outgoing the 
star opens and subsequently varies about an average value. 
To determine this magnetic flux we calculated the depen- 
dence of the maximum value of the magnetic flux function 
&(r, 6,t) = r sin 9A^, = maxg^ro, 0, t) on time at the 
sphere of radius ro = 30. The magnetic field lines along 
which ^ > arc evidently closing inside the sphere of 
radius r . The value \&o determines the location of the last 
closed field line at this distance, in particular the radius at 
which this field line crosses equatorial plane r eq = p/^ {t). 
That is, r e q shows the radius of the closed magnetosphere. 
The right-hand panel of Figure 7 shows that this radius 
r eq {t) varies with time in a fashion similar to rd(t); how- 
ever, the amplitude is much smaller. 

An important aspect of the interaction between a 
rapidly rotating magnetized star and an accretion disk is 
the transport of angular momentum. The equation for the 
conservation of angular momentum can be obtained from 
(16) by multiplying it by rsin#, 



d(pv^r sinff) 1 d 



dt 

1 o> 
r sin 9 86 



dr 



r sin 9[ T, 



v t pr sinfl— - 
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sin 9 



rsm( 



Ooj 

Te<f> - f t p sin 0— 



= 



where T r(j> = pv r v^ - BrB^/An and T g<ji = pvgv^ ~ 
BgB^/An are the components of "non- viscous" part of the 
stress tensor. 

The angular momentum conservation equation evidently 
has the form of a continuity equation, 



d(pi) 
dt 



V-L = , 



where I = v^r sin 9 is the specific angular momentum and 
the angular momentum flux density is 



L = r sin I pv^v v 



B^Hp 
4ir 



v t pr sin 9X7 oj 



4. RESULTS 

We performed simulations of accretion to a star in the 
"propeller" regime for a wide range of parameters, \i, £1*, 



The first term on the right-hand-side L gives the trans- 
port of angular momentum by the matter; the second term 
is magnetic field contribution; and the third term is the 
transport due to viscous stress. 
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Fig. 6. — The left-hand panel shows the inner part of the simulation region with the logarithm of density as a background 
and the magnetic field lines (solid lines) and the grid. The right-hand panel shows the full simulation region. 



Figure 8 (a,b,e-g) shows the angular momentum fluxes 
at t = 2300. The streamlines show the direction of the 
fluxes, while the background shows the magnitude of the 
flux. Figures 8 (a, e) show the distribution of the angular 
momentum fluxes carried by the magnetic field at large 
and small scales. 

Figures 8(b, f) show the corresponding fluxes carried 
by the matter. From Figure 8e it is clear that the an- 
gular momentum outflow from the star is due mainly to 
the twisting of the magnetic field. Note that a larger part 
of angular momentum is carried along the open field lines. 
Other part of angular momentum which is along the closed 
field lines is transported from the star to the disk. Angular 
momentum carried by matter is largest in the disk and in 
the conical outflow in the vicinity of the neutral layer (Fig. 
8b). Here magnetic pressure is not very high, while pres- 
sure and density of matter play larger role. Bottom panels 
show energy fluxes carried by the magnetic field (Fig. 8c) 
and by matter (Fig. 8d). One can see from Figure 8c that 
significant energy is carried from the star by a Poynting 
flux. 

Figure 8g shows the angular momentum flux carried by 
the viscous stress. This stress is large at the boundary 
between the disk and magnctosphcre in the vicinity of the 
equatorial plane. This is the place of significant gradients 
of the angular velocity of matter, as one can see from the 
Figure 9 (shown at the same time t — 2300). This figure 
also shows the distribution of density, entropy (logarithm 
scale at the right), and radial velocity v r in the equatorial 
plane. It is evident that the disk has a sharp boundary 
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Fig. 7. — The left-hand panel shows the variations of the 
variations of the radius where a given magnetic field line cro 



both in density and in entropy. 

Next, we evaluate the directions where most of the mat- 
ter and angular momentum outflows. For this we calcu- 
lated the distribution of the fluxes of matter and angular 
momentum as a function of 9 (the colatitude) as a func- 
tion of time at a given distance tq = 20. Figure 10 shows 
a plot of these distributions as a function of time. One 
can see that the matter flux to the wind is mainly in the 
range 40° < 9 < 70°. Thus, most of the matter flows into 
a wide-angle hollow-cone with a large inclination angle to 
the z— axis (Figure 10a). There is also a flux of matter 
towards the star which is associated with the equatorial 
region of the disk. The angular momentum flux carried 
by matter has distribution similar to the distribution of 
the matter flux (see Figure 10b). The angular momentum 
flux carried by the magnetic field (Figure 10c) is concen- 
trated in two stripes: 0° < 9 < 30° (near the axis), and 
60° < 9 < 80° (closer to the disk) . The larger flux is asso- 
ciated with the angular momentum outflow in collimated 
jet, along the field lines which start at the surface of the 
star. The other stripe, at larger angles 9 is connected with 
the disk wind. 

The wide-angle, hollow-cone disk wind is similar to the 
stationary magnetocentrifugally driven winds predicted by 
Blandford and Payne (1982) and Lovelace, Berk, and Con- 
topoulos (1991) and first obtained and analyzed in MHD 
simulations by Ustyugova et al. (1999). For such winds 
the angle of the poloidal field to the z— axis at the sur- 
face of the disk O = ta,n~ 1 (\B r /B z \) is predicted to be 
^> 30°. On the other hand a magnetically dominated (or 
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inner radius of the disk. The right-hand panel shows the 
ses the equatorial plane. 
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Fig. 8. — The top panels show the angular momentum flux carried by the field (a) and by the matter (b). The bottom 
panels show the energy fluxes carried by the field (c) and by the matter (d). The right-hand panels show the angular 
momentum fluxes close to the star carried by the field (e), by the matter (f), and by the viscous stress (g). The background 
color indicates the value of the flux while the streamlines indicate the direction of the flux. 




Fig. 9. — Radial distribution of different variables in the equatorial plane, the density p, the specific entropy 5, the radial 
velocity v r , the angular rotation rate u> = v^/r, and the Keplerian angular velocity lok- The scale for the entropy is on 
the right-hand side. 



Poynting) outflow from a disk or star can have the angle 
9 less than 30° (Ustyugova et al. 2000; Lovelace et al. 
2002). The narrow-angle, hollow-cone jet along the open 
field lines of the star has much smaller values of O and it 
is magnetically dominated. 

Figure 11 shows the distribution of different physical 
quantities in the computational region (top panels) and 
on a smaller scale (bottom panels). Figures 11 (a, d) show 
distribution of the matter flux pv p , and streamlines of mat- 
ter flow at time t = 2300. One can see that matter flows 
with the disk, then turns near the magnetosphere of the 
star and flows outward in a disk wind. Some matter goes 
around the magnetosphere and flows to the jet along the 
field lines starting at the stellar surface. The white line 



shows the escape velocity v esc — ^2GM jr. Above this 
line v > v esc . The color background of Figures 11 (b, e) 
shows the plasma (3 = 8irp/B 2 parameter. The dashed 
white line corresponds to j3 = 1. Figures 11 (c, f) show 
isocontours of the r sin OB^ which is the poloidal current 
flow through a horizontal disk < r sin 6 centered on the 
z— axis. 

Figure 12 shows plots of different variables as a func- 
tion of distance along a magnetic field line which starts 
from the disk in the vicinity of the neutral field line where 
the matter flux to the wind is the largest. The solid line 
shows the projection of the poloidal velocity onto this field 
line, the dashed line shows the projection of the velocity 
component normal to the magnetic field line. One can 
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Fig. 10. — Angular (6) distributions of different fluxes on the surface of a sphere of radius tq = 20 as a function of time: 
(a) shows the matter flux; and (b) and (c) show the radial angular momentum flux carried by the matter and by the 
magnetic field, respectively. 
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Fig. 11. — In panel (a) the background shows the matter flux and the lines are the streamlines of the matter flow. The 
white line is the v = v esc line. In panel (b) the background shows the plasma parameter (3 = 8np/B 2 . The streamlines 
are the magnetic field lines, and the white dashed line shows the /3 = 1 line. In panel (c) the background shows rsin^i?^, 
which for steady conditions is the total poloidal current through a disk < rsin# center on the z— axis. The lines are the 
magnetic field lines. The white dashed line shows the neutral line where B p — 0. The bottom panels (d - f) show the 
same as (a - c) but in the vicinity of the star. 



see that the poloidal velocity is nearly parallel with the 
poloidal field, which is in accord with the theory of sta- 
tionary axisymmetric flows of ideally conducting plasma 
(see, e.g., Ustyugova et al. 1999). The sound speed c s and 
escape velocity v esc are also shown. One can see that the 



flow becomes supersonic and that the flow velocity exceeds 
the escape velocity v esc . Earlier, Ustyugova et al. (1999) 
obtained dependences similar to those of Figure 12 for the 
case of an ordered magnetic field threading a disk around 
a non-magnetized star. 
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Fig. 12. — Velocities along a field line which starts from the disk near the neutral line of the magnetic field (where the 
matter flux to the wind is the largest). v\\ is the velocity component parallel to the field line and v± is the component 
perpendicular to the field line. v esc is escape velocity and c s is a sound speed. 



Interaction of the magnetized star with the accretion 
disk can lead to the outflow and escape of a significant 
fraction of the matter incoming in the disk. We find in 
general that there is a jet and a disk wind giving outflows 
of energy, angular momentum, and matter from the sys- 
tem. The jet is identified as the flow in the region within 
the neutral line of the poloidal magnetic field. The mag- 
netic flux in the jet is the flux emanating from the star. 
The disk wind is the flow outside of the neutral line of the 
poloidal magnetic field. For both outflows we take into 
account only the matter which has a velocity higher than 
the escape velocity v esc . 

Thus, the simulation region includes three flow regions, 
the disk, the jet, and the wind. These regions do not 
constitute the entire simulation volume because part of it 
is filled with relatively rarefied matter of magnetosphcrc 
and corona. Accreting matter which is not accelerated to 
escape speed goes into the magnetosphere or the corona. 
These regions do not have a significant role in the mat- 
ter and angular momentum flow. The star has an impor- 
tant role in absorbing mass and losing angular momentum. 
The star influences the plasma dynamics through the inner 
boundary of the simulation region. However the influence 
of the plasma on the star is noticeable on only an extremely 
long time scale. Figure 13 shows the interactions between 
the different elements of the system. Evidently, matter 
which comes through an accretion disk flows in part to the 
disk wind, in part to the jet, and in part accretes to the 
star. There is also angular momentum transport between 
different elements. The star represents a large reservoir 
of angular momentum which is transported to the disk 
through the closed field lines and to the jet through the 
open field lines. The positive direction of the fluxes is 
shown by arrows. 

In all cases the fluxes of matter and angular momentum 
oscillate strongly in time, but always vary around an aver- 
age value. See for example Figure 10. For this reason we 
calculated the time-averaged fluxes and investigated their 
dependences on the parameters of the model. We took as 
a base the mentioned reference case with dimensionless pa- 
rameters jU = 10, fi* = 1, a v = 0.2, and ay = 0.2. A series 
of simulations was then done to obtain the dependences 
of the time-averaged fluxes on different parameters. Only 
one parameter was varied in each series. Figures 14-16 



show the derived dependences. 

Figure 14 (a-d) show the dependences of the matter 
fluxes on our main parameters. Here, M<j is the mass 
accretion rate in the disk at the radius r$ = 30, M s is 
the accretion rate to the star, M w is the matter flux to 
the wind, Mj is the matter flux to the jet. As we men- 
tioned, we distinguish between "strong" (with outflows) 
and "weak" (no outflows) propellers. Results of model- 
ing of the "strong" propellers are marked with the filled 
symbols, while the "weak" propeller results are indicated 
by open symbols. Different symbols show different depen- 
dences: Mj (squares), M w (circles), M4 (triangles), and 
M s (gradient signs). The solid lines show approximations 
of these dependences. Most of dependences could be ap- 
proximated by a simple power law and have a clear sense. 
For example, the larger the magnetic moment /i and the 
larger the angular velocity of the star f2* , the larger matter 
fluxes to the wind and jet and the smaller the accretion 
rate to the star. At smaller fi and ejection of matter 
to the wind and jet becomes less efficient and we enter the 
"weak" propeller regime with very weak or no outflows 
(see Fig. 14 a,b). Note, that all "weak" propellers arc on 
the left-hand side of the plots. 

It is interesting to look at dependences on viscosity, a v 
(see Fig. 14c). We see that matter fluxes to the wind/jet 
strongly increase with increase of viscosity, while the ac- 
cretion rate to the star decreases. We conclude that the 
viscous stress contributes to launching the outflows. At 
sufficiently small viscosity, a v <; 0.1, outflows are absent. 
The role of the viscosity is twofold. From one side, at 
larger viscosity, the "friction" between the inner regions 
of the disk and magnetosphere is larger. From other side, 
the radial velocity in the disk is proportional to a v , so that 
viscosity increases the matter flux per unit area and the 
inner region of the disk come closer to the inner regions 
of the fast rotating magnetosphere. Both factors lead to 
enhanced outflows, but it is difficult to separate the two 
factors. In RUKL05 test simulations were done for a case 
with a v = 0.1, and with the density in the disk doubled. 
This of course increased the accretion rate but it also en- 
hanced the outflows. 

Figure 14d shows the dependence on the magnetic dif- 
fusivity. The dependence of the matter fluxes to the 
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Fig. 13.— Scheme of interaction between different elements of the system. 
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Fig. 14.— Variation of the matter fluxes as a function of different parameters. 



wind/jet on ad is more complicated than other depen- 
dences. Namely, for ad <J 0.2, the fluxes M w and Mj de- 
crease with ad, because the mixing of the disk matter with 
the magnetic field of the magnetosphere become less and 
less efficient, and correspondingly the angular momentum 
transport from the star to the disk matter decreases. On 
the other hand, for ad ;> 0.2, the diffusivity become too 
high, and the inner regions of the disk and magnetosphere 
are not coupled sufficiently to transport angular momen- 
tum. This explains the parabolic approximations of these 
dependences. We took into account all dependences shown 
in Figures 14 (a-d), and approximated them with analytic 
functions. We approximated the dependences on ad with 
power laws to the left and to the right of the dividing value 
ad ~ 0.2. In Figures 14 (c-d), we show the dependences 
around a v ~ 0.2 and ad ~ 0.2. We have made additional 



runs not shown in the plot at much smaller a v and ad, 
down to 0.01. Also we have made runs for much larger 
values, up to 0.6. We show below dependences which in- 
corporated all these runs, separately, for ad > 0.2 and 
ad < 0.2. The matter fluxes to the wind M w , to the jet 
Mj, through the disk Mj, and to the star M s were approx- 
imated as: 
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*-~-(4r(sr<g) M . 

The matter flux to the disk wind increases with /i al- 
most linearly and strongly increases with f2*, while the 
accretion rate to the star M s strongly decreases. We note 
that dependences on /i and f2» have a threshold charac- 
ter. Namely, the matter fluxes to the wind and to the jet 
strongly decrease for fi <; 7 and fi* <; 0.8, where the sys- 
tem enters the "weak" propeller regime (see open symbols 
in Figures 14 a,b). One can see from the Figure 14b that 
the dependence of M w and Mj on the angular velocity of 
the star fi* also has a threshold behavior for fi* <; 0.8. 
These threshold dependences reflect the qualitative dif- 
ference between "strong" propellers (with outflows) and 
"weak" propellers (with no outflows). 

As a measure of the efficiency of the propeller interac- 
tion with the disk, we take the ratio of the outflow rate 
to the disk wind and jet to the accretion rate through the 
disk, S = Mj W /Md- For the larger values of S, less matter 
reaches the surface of the star. Correspondingly, a larger 
part of the disk matter acquires additional angular mo- 
mentum and flows to the jet and/or wind. Figures 15 a-d 
show the dependences of the value S on the parameters of 
the model. The efficiency of the propeller depends weakly 
on the magnetic moment of the star. The larger the angu- 
lar velocity of the star and the larger the a-coefficient of 
viscosity, the higher is the efficiency of the propeller. The 
dependence of the efficiency on the a-coefficient of mag- 
netic diffusivity has a maximum near ad ~ 0.2. We found 
the following dependences, 
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Another measure of the efficiency of the propeller is the 
ratio of the total matter flux going to the wind /jet to the 
matter flux to the star, 
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The ratio of outflow to accretion increases strongly with // 
and fi*. 

Next, we discuss angular momentum transport from the 
star and its dependence on main parameters of the model. 
The outflow of angular momentum from the star causes it 
to spin down. Part of the angular momentum is carried by 
the matter, another part by the tension of the magnetic 
field lines, and a further part by the viscous stress. 

The flux of angular momentum from the star was cal- 
culated as 
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where the integration was over the surface of the inner 
boundary of the simulation region which is close to the 



surface of the star and where dS is the surface element di- 
rected outward to the region. Simulations show that close 
to the star angular momentum flux is determined mainly 
by the magnetic stress. 

Part of the flux is associated with the closed field lines. 
These are field lines connecting a star and the disk. Other 
part is associated with the open field lines, which connect 
a star with the low-density corona. We calculated a total 
flux N 8 and also a flux associated with the closed field lines, 
N s d, that is angular momentum flux transported from the 
star to the disk. 

Figure 16 shows fluxes N s (gradient signs) and N s d (cir- 
cles). We found the following dependences on different 
parameters: 
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Figure 16 (a, b) shows that both fluxes increase with fi 
and fi*. They increase with a similar rate so that the lines 
in Figures 16 (a,b) are almost parallel. Figures 16 (c,d) 
show that the fluxes depend only weakly on the a coeffi- 
cients of viscosity and diffusivity. They have maximum at 
a d ~ 0.2. Note that the angular momentum flux carried 
by the open field lines of the jet N s — N s d is similar to that 
carried by the closed field lines. That is, a star in the pro- 
peller regime in our range of parameters spins-down due 
to both: open and closed field lines. 

Now, we can estimate the time-scale of spin-down for 
CTTSs and accreting neutron stars for parameters given in 
§3.4. For CTTSs, the loss of the angular momentum of the 
star in our main case is N s — 3.1Lo — 1.06 x 10 37 gcm 2 /s 2 . 
The star's angular velocity is fi = 2tt/P w 7 x 10~ 5 s _1 , its 
angular momentum is J = kMr 2 Vt = 2.2 x 10 51 fc gcm 2 /s, 
where k < 1. Taking k = 0.4, the spin-down time-scale is 
r = J/N s w 2.7 x 10 6 years. Note, that this time-scale was 
calculated for a magnetic field £?* = 10 3 G and a relatively 
low accretion rate, M w lO _8 M0/yr. The time-scale de- 
creases with magnetic field of the star as <~ B~ 1A and will 
be r s=s 8.1 x 10 5 years for B„ ~ 3 x 10 3 G. Time-scale also 
decreases with matter flux. Thus, if CTTSs have a strong 
magnetic field in the past, then they were at the propeller 
regime and were spun-down in < 10 6 years. This time- 
scale is shorter than typical life-time of the CTTSs, which 
is ~ 10 7 — 10 8 years. We conclude that propeller mecha- 
nism may be responsible for fast spinning-down of CTTSs 
to presently observed slow rotation in the early stages of 
their evolution. 

For an accreting neutron star, N s = 3.1L = 3.87 x 
10 33 gcm 2 /s 2 , fl = 2tt/P w 5 x 10 3 s" 1 . The star's angular 
momentum is J = kMr 2 il = 1.5 x 10 49 fc gcm 2 /s. The 
time-scale of spin-down is r = J/N s w 2.5 x 10 7 years. 
This time-scale is shorter than expected life-time of mil- 
lisecond pulsars in the accreting stage. However, millisec- 
ond pulsars have a different history of evolution compared 
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Fig. 15. — Efficiency of the propeller as a function of different parameters as discussed in the text. 
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Fig. 16.— Variation of the angular momentum fluxes with parameters of the system . 



to CTTSs, and they expected to be spined up during ac- 
creting stage (Bisnovatyi-Kogan & Romberg 1974). How- 
ever, the pulsar mechanism may be responsible for fast 
spinning-down of the pulsar during periods of decreased 
accretion rate. 

5. CONCLUSIONS 

We performed a large number of axisymmetric MHD 
simulations and did extensive analysis of the interaction of 



a rapidly rotating magnetized star with an accretion disk 
in the propeller regime. The results can be divided to two 
classes, "strong" propellers and "weak" propellers. The 
"strong" propeller is characterized by the intense, wide- 
angle, conical matter outflow in the vicinity of the neutral 
line of the poloidal magnetic field and formation of well col- 
limated, magnetically dominated jets along the rotational 
axis. In the "weak" propeller regime there are no appre- 
ciable matter outflows from the disk and the collimated 



18 



flow along the axis is much weaker. 

In this work the four main parameters of the model were 
varied: the magnetic moment of the star, its angular ve- 
locity, and a— coefficients of the kinematic viscosity and 
magnetic diffusivity. We calculated time-averaged total 
fluxes of matter, energy, and angular momentum and in- 
vestigated their dependences on the main parameters of 
the model. The average characteristics were calculated 
during 1000 — 2500 periods of rotation of the inner disk. 
The derived dependences were approximated in most cases 
by power law functions. But some dependences showed a 
threshold character with the flux becoming negligible for 
a range of the parameters. 

In all cases the interaction of the magnetized star with 
the accretion disk showed strongly non-stationary, quasi- 
periodic character. For example, in the reference case de- 
scribed in the paper, the typical "period" of oscillations is 



several tens of our time-units. All fluxes of matter, energy, 
and angular momentum to the wind, jet, and star oscillate 
strongly with time. 

Our simulation results are applied to classical T Tauri 
stars and to millisecond pulsars. We conclude that pro- 
peller mechanism may be responsible for fast spinning- 
down of CTTSs to presently observed slow rotation dur- 
ing the first <~ 10 6 years of their evolution. In the case 
of accreting, magnetized millisecond pulsars, the propeller 
mechanism may be responsible for relatively fast spin- 
down during periods of lower accretion. 
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